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Abstract 

In this paper we continue our studies of the phase space geometry and dynamics associated with 
index k saddles {k > 1) of the potential energy surface. Using Poincare-Birkhoff normal form 
theory, we give an explicit formula for a "dividing surface" in phase space, i.e. a co-dimension 
one surface (within the energy shell) through which all trajectories that "cross" the region of the 
index k saddle must pass. With a generic non-resonance assumption, the normal form provides 
k (approximate) integrals that describe the saddle dynamics in a neighborhood of the index k 
saddle. These integrals provide a symbolic description of all trajectories that pass through a 
neighborhood of the saddle. We give a paramctrization of the dividing surface which is used as 
the basis for a numerical method to sample the dividing surface. Our techniques are applied to 
isomerization dynamics on a potential energy surface having 4 minima; two symmetry related pairs 
of minima are connected by low energy index one saddles, with the pairs themselves connected via 
higher energy index one saddles and an index two saddle at the origin. We compute and sample 
the dividing surface and show that our approach enables us to distinguish between concerted 
crossing ("hilltop crossing") isomerizing trajectories and those trajectories that are not concerted 
crossing (potentially sequentially isomerizing trajectories). We then consider the effect of additional 
"bath modes" on the dynamics, which is a four degree-of-freedom system. For this system we 
show that tlic normal form and dividing surface can be realized and sampled and that, using the 
approximate integrals of motion and our symbolic description of trajectories, we are able to choose 
initial conditions corresponding to concerted crossing isomerizing trajectories and (potentially) 
sequentially isomerizing trajectories. 

PACS numbers: 05.45.-a, 45.10.Na, 82.20.Db, 82.30.Qt 
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I. INTRODUCTION 



Transition state theory has long been, and continues to be, a cornerstone of the theory of 
chemical reaction rate^^. A large body of recent research has shown that index one saddled 
of the potential energy surfac^^ give rise to a variety of geometrical structures in phase 
space, enabling the realization of Wigner's vision of a transition state theory constructed in 
phase spac^^^. 

Following these studies, attention has naturally focussed on phase space structures asso- 
ciated with saddles of index greater than one, and their possible dynamical significanc^^^^. 
In previous work we have described the phase space structures and their influence on trans- 
port in phase space associated with index two saddles of the potential energy surface for n 
degree-of-freedom (DoF) deterministic, time-independent Hamiltonian system^. (The case 
of higher index saddles has also been investigated by Haller et ahl^SEIl gee also refs [3514371 ) 

The phase space manifestation of an index one saddle of the potential energy surface in 
an n DoF system is an equilibrium point of the associated Hamilton's equations of saddle- 
center-. . .-center stability type. This means that the matrix associated with the linearization 
of Hamilton's equations about the equilibrium point has one pair of real eigenvalues of equal 
magnitude, but opposite in sign (±A) and n — 1 pairs of purely imaginary eigenvalues, iiwj, 
j = 2, . . . , n. 

The phase space manifestation of an index k saddle is an equilibrium point of saddle 
stability type: 

saddle x . . . x saddle x center x . . . x center . (1.1) 

k times n—k times 

The matrix associated with the linearization of Hamilton's equations about the equilibrium 
point then has k pairs of real eigenvalues of equal magnitude, but opposite in sign (±Aj, 
i = 1, . . . , /c) and n — k pairs of purely imaginary eigenvalues, ±ioOj, j = /c + 1, . . . , tW^. 
Informally, an index k = 2 saddle on a potential surface corresponds to a maximum or 
"hilltop" in the potentiaPSEH. 

Although it has been argued on the basis of the Murrell-Laidler theorenP^Sl \}igx critical 
points of index 2 and higher are of no direct chemical significanc^^'^, many instances of 
index 2 (and higher) saddles of chemical significance have been identifiecP^. For example, 
Heidrich and Quapp discuss the case of face protonated aromatic compounds, where high 
energy saddle points of index two prevent proton transfer across the aromatic ring, so that 
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proton shifts must occur at the ring periphery^. Index two saddles are found on potential 
surfaces located between pairs of minima and index one saddles, as in the case of internal 
rotation/inversion in the H2BNH2 molecul^^or in ureaP', or connected to index one saddles 
connecting four symmetry related minima, as for isomerization pathways in 

Saddles with index > 1 might well play a significant role in determining system properties 
and dynamics for low enough potential barrierpSMI or at high enough energie^. The role 
of higher index saddles in determining the behavior of supercooled liquids and glasseP^^ is 
a topic of continued interest, as is the general relation between configuration space topology 
(distribution of saddles) and phase transitions^. 

Several examples of non-MEP (minimum energy path) reactionS ^^ * ^^^ and "roaming" 
mechanism^S^'^ have been identified in recent years; the dynamics of these reactions is 
not mediated by a single conventional transition state associated with an index one saddle. 
Higher index saddles can also become mechanistically important for structural transforma- 
tions of atomic cluster^ when the range of the pairwise potential is reducecP^. We note 
for example the work of Shida on the importance of high index saddles in the isomerization 
dynamics of Ki-j clusterd^. 

The role of index two saddles in the (classical) ionization dynamics of the Helium atom 
in an external electric fielcP^^ has recently been studied from a phase space perspective by 
Haller et al.!™. 

In previous work we have discussed phase space structures and their influence on phase 
space transport in some detail for the case of an index two saddle of the potential energy 
surface corresponding to an equilibrium point of saddle-saddle-center-. . .-center stability 
typd^. In this paper we extend our analysis in several respects. 

One motivation for the work presented here is the possibility of developing a dynamically 
based characterization of concerted and sequential reaction pathwayd^SES! in phase space. 
Consider isomerization dynamics on the model potential shown in Fig. [l] (This potential 
is discussed in more detail in Sec. |Vj) The potential has 4 minima; two symmetry related 
pairs of minima are connected by low energy index one saddles, while the pairs themselves 
are connected via higher energy index one saddles. An index two saddle (hilltop, denoted 
II) is located at the origin. In Figure [T| we indicate schematically two possible pathways 
between the lower left hand well (designated ( ); the symbolic code is discussed further in 



Sec. HI) and the upper right hand well, (++) (see also Fig. 7 of ref. [THI). At energies above 
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that corresponding to the index two saddle, there are, quahtatively speaking, two possible 
isomerization routes: a sequential path, shown in Fig. [T^, in which a trajectory passes 
through an intermediate well (either (H — ) or ( — h)) in the course of the isomerization, and 
a concerted route, shown in Fig. in which the trajectory effectively passes directly from 
the reactant well to product well without entering a well corresponding to an 'intermediate' 
species. 

It is natural to ask whether the above qualitative distinction between sequential and 
concerted reaction pathwayd^SBSl, made on the basis of considerations of the nature of iso- 
merizing trajectories in configuration space, can be given a rigorous dynamical formulation 
in phase space. In the present work we provide such a formulation, based on the properties 
of the normal form in the vicinity of the index two saddle point. In particular, we show 
that it is possible to define phase space dividing surfaces for higher index saddles which 
generalize the now-familiar dividing surfaces defined for index one saddlePHsHto the case of 
higher indices, and that phase points can be chosen on such dividing surfaces with prescribed 
dynamical character (e.g., concerted crossing trajectories). 

The structure of this paper is as follows. In Sec. [TT] we review and extend our previous 
analysi^ of the phase space structure in the vicinity of a (non-resonant) index k saddle in 
terms of the normal form. Particular emphasis is given to discussion of the dynamical sig- 
nificance of the values of the associated action integrals, and to the symbolic representation 
of the qualitatively distinct classes of trajectory behavior in the vicinity of the saddle-saddle 
equilibrium. For generic non-resonance conditions on the eigenvalues of the matrix associ- 
ated with the linearization of Hamilton's equations about the equilibrium point, the normal 
form Hamiltonian is integrable. Integrability provides all of the advantages that separability 
provides for quadratic Hamiltonians: the saddle dynamics can be described separately and 
the integrals associated with the saddle DoFs can be used to characterize completely the 
geometry of trajectories passing through a neighborhood of the equilibrium point. As for the 
case of index one saddles, normally hyperbolic invariant manifolds (NHIMs)?^'^ associated 
with index two saddles are an important phase space structure and we give a brief discussion 
of their existence and the role they play in phase space transport in the vicinity of index k 
saddles, following our work irP^. 



In Sec. Ill we introduce the concept of concerted crossing trajectories. The concerted 
crossing trajectories are defined in terms of the symbolic code introduced previouslj^, and 
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realize the inuitive notion of direct, concerted passage between wells via the hilltop region. 
Those trajectories that are not concerted crossing are potentially, but not necessarily, se- 
quentially isomerizing trajectories. 



Sec. IV defines the dividing surface (DS) for index k saddles. This is a key aspect of 
the present work; we show that it is possible to define a codimension one surface in phase 
space through which all concerted crossing trajectories must pass, and which is everywhere 
transverse to the fiow. We also introduce a parametrization of the dividing surface. This 
parametrization together with the use of the normal form enables us to sample the DS and 
select phase points on trajectories having specified dynamical character. 

In Section |V] we present a numerical study of the phase space structure in the vicinity of 
an index two saddle in the context of a problem of chemical dynamics, namely, isomerization 
in a multiwell potential. For isomerization on this model potential energy surface, which 
has multiple (four) symmetry equivalent minima, analysis of the phase space structure in 
the vicinity of the index two saddle enables a rigorous distinction to be made between con- 
certed crossing ("hilltop crossing") isomerizing trajectories and those trajectories that are 
not concerted (potentially sequentially isomerizing)^^'^. Our normal form based procedure 
for sampling the DS enables us to determine phase points lying on concerted crossing trajec- 
tories. Numerical propagation forwards and backwards in time shows that such trajectories 
do indeed have the properties intuitively associated with the concerted isomerizing pathway. 
Sec. |Vl] concludes. 
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II. THE POINCARE-BIRKHOFF NORMAL FORM IN A PHASE SPACE 
NEIGHBORHOOD OF AN INDEX k SADDLE 



We begin by considering the normal form in the neighborhood of an equihbrium point 
of saddle- . . .-saddle-center-. . .-center stability type, where there are k saddle degrees-of- 
freedom (DoF) and n — k center degrees-of-freedom (eq. )). We assume the usual non- 



resonance condition on the eigenvalues of the matrix associated with the linearization of the 
Hamiltonian vector field about the equilibrium point (see, e.g.^. In particular, we assume 
that the purely imaginary eigenvalues satisfy the non-resonance condition k^^iUk+i + . . . + 
kn^n 7^ for any {n — A;)-vector of integers (fcfc+i, • • • ,kn) with not all the ki = (that 
is, {kk+i, . . . , kn) G Z"^^ — {0}) and the real eigenvalues satisfy the (independent) non- 
resonance condition kiXi + . . . + kkXk for any A;-vector of integers (fci, . . . , kk) with not 
all the ki = 0, i = 1, . . . , k. In this case the normal form transformation transforms the 
Hamiltonian to an even order polynomial in the variables 

li = QiPi, i = 1,. . . ,k, (2.1a) 
I, = ^{q]+p]),J = k + l,...,n. (2.1b) 

In other words, we can express the normal form Hamiltonian as: 

i7(/i,J2,J3,...,/„), (2.2) 

with associated Hamilton's equations: 

OH dHdIi , 
gi = 7^ = TT^T^ = Ai^i, (2.3a) 
dpi dli dpi 

dH dHdIi ^ , , , , 

Pi = = = -AiPi, 1 = 1, ...,k, (2.3b) 
dqi dli dqi 

Pj = -— = - — — = -VLjq^, J = k + l,...,n, (2.3d) 



dqj dij dqj 
where we have defined the frequencies 



f)TJ 

A,(X)^— (X), t = l,...,k, (2.4a) 
dli 

n,{X) = —{X), j = A; + l,...,n, (2.4b) 
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where X = {Ii, I2, 13, . . . , In) and it can be verified by a direct calculation that the Ij, 



j = 1, . . . ,n, are integrals of the motion for (2.3). 



The integrability of the normal form equations (2.3) provides us with a very straightfor- 
ward way of characterizing the dynamics in a neighborhood of the index k saddle. The coor- 
dinates {qi,Pi), i = I, . . . , k, describe saddle-type or 'reaction dynamics' , and the dynamics 
in the k saddle planes can be further characterized by the saddle integrals, li, i = 1, . . . , k. 
The coordinates {qi,Pi), i = k + 1, . . . ,n, describe bounded motions (center-type dynamics) 
or 'bath modes', which can be further characterized by integrals li, i = k + 1, . . . ,n. 

We now introduce a canonical transformation of the saddle variable^. Passage of tra- 
jectories over the saddle is more naturally described in terms of these new coordinates, and 
computation of codimension one dividing surfaces for index k saddles is also facilitated. The 
transformation is given by: 

(li = —7= {Qi + Pi) , (2-5a) 



V2 
1 
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Pi = A= {pi - Qi) , i = l,...,k, (2.5b) 



with inverse 



= ^ {Qi - Pi) , (2.6a) 
Pi = -^{Pi + Qi) , i = l,...k. (2.6b) 



The transformation of variables given by eq. (2.6), where i = 1, . . . , fc, and the remaining 
variables are transformed by the identity transformation, is canonical. The variables g^, 
i = 1, . . . ,k, are naturally identified with physical configuration space coordinates in the 



vicinity of the saddle. The Hamiltonian is given by eq. (2.2), with action variables 

li = QiPi = ^{p^i -Qi) ,i = l,...,k, (2.7a) 
Ij = \{q]+p',),3 = k + l,...,n, (2.7b) 
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and Hamilton's equations then take the following form: 

, OH dHdIi ^ _ , , 

OH dHdIi ^ _ , , , , 

~~dq- ^ ^~dT~dq-^ ~ « = l,...,fc, (2.8b) 

dH dHdlj ^ 
p, = = = -^-^.^ J = k + l,...,n, (2.8d) 

A. A Normally Hyperbolic Invariant Manifold 

As noted previousl}^, for n > /c the 2n — 2k — 1 dimensional surface: 

M = {qi=Pi = --- = qk=Pk = 0, H{0, . . . , 0, h+i, ...,/„) = E > 0} , (2.9) 

is a normally hyperbolic invariant manifold (NHIM) in the energy surface H{Ii, . . . , /„) = 
E > 0. Moreover, this NHIM has 2n — k — 1 dimensional stable and unstable manifolds 
(within the fixed energy surface). Note that these stable and unstable manifolds are codi- 
mension one in the energy surface only for k = 1, i.e., index one saddles. Nevertheless, in 
the construction of dividing surfaces for index k saddles we will see that the NHIM (but not 
its stable and unstable manifolds) plays a similar role for all k. 

The case n = k, i.e. a n degree-of-freedom system with an index n saddle requires special 
consideration. In this case the NHIM corresponds to an equilibrium point and exists only 
on the E = energy surface. Thus, to get a non-trivial NHIM existing for i^^ > we must 
have n > k. 



B. Accuracy of the Normal Form 

Here we briefly address some of the issues associated with the accuracy of the normal 
form. The transformation to normal form is implemented by an algorithm that operates 
in an iterative fashion by simplifying terms in the Taylor expansion about the equilibrium 
point order by order, i.e., the order M terms are normalized, then the order M + 1 terms 
are normalized, etc!^. The algorithm is such that normalization at order M does not affect 
any of the lower order terms (which have already been normalized). The point here is that 
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although the algorithm can be carried out to arbitrarily high order, in practice we must stop 
the normalization (i.e., truncate the Hamiltonian) at some order M, after which we make a 
restriction to some neighborhood of the saddle in which the resulting computations achieve 
some desired accuracy. It is therefore necessary to determine the accuracy of the normal form 
as a power series expansion truncated at order M in a neighborhood of the equilibrium point 
by comparing the dynamics associated with the normal form to the dynamics of the original 
system. There are several independent tests that can be carried out to verify accuracy of 
the normal form, such as the following: 

• Examine the extent to which integrals associated with the normal form are conserved 
along trajectories of the full Hamiltonian (the integrals will be constant on trajectories 
of the normal form Hamiltonian). 

• Check invariance of specific invariant manifolds (e.g., the NHIM, its stable and unstable 
manifolds, the energy surface) under dynamics determined by the full Hamiltonian. 

Both of these tests require us to use the transformation between the original coordi- 
nates and the normal form coordinates. Software for computing the normal form as well 
as the transformation (and inverse transformation) between the original coordinates and 
the normal form coordinates can be downloaded from |http : / /lacms . maths . bris . ac . uk/ 



publications/software/index.html Specific examples where the accuracy of the normal 
form and its relation to M, the fixed neighborhood of the saddle, and the constancy of 
integrals of the truncated normal form on trajectories of the full Hamiltonian can be found 
in refs [T3p5] - [T71l78l A general discussion of accuracy of the normal form can be found in ref. 
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III. CROSSING AND CONCERTED CROSSING TRAJECTORIES 



In this section we introduce the notion of crossing trajectories. Crossing trajectories pass 
through a neighborhood of the index k saddle in such a way that all the saddle coordinates % 
change sign, and constitute a dynamically well-defined subset of trajectories. Our treatment 
of crossing trajectories provides an excellent illustration of the power of the normal form 
(NF) and the utility of the integrals of the motion in analyzing dynamics in the vicinity of 
the saddle. 

The definition of crossing trajectories given here is purely local, in that it relies on the 
values of the integrals computed using the NF in the vicnity of the saddle. In our discussion of 
the isomerization dynamics for a model multi-well system given below in Sec. |V| we establish 
a connection between the local crossing property and the more global mechanistic notion of 
concerted (as opposed to sequential) isomerization trajectories; in this context, the crossing 
trajectories defined here are usefully described as concerted crossing (CC) trajectories. 

In the normal form coordinate system the 'crossing' of a saddle is analyzed by considering 
only the saddle degrees-of-freedom, {qi,pi), i = l,...,k, since the remaining coordinates 
remain bounded. Crossing occurs when all coordinates qt, i = 1, . . . , k, change sign as the 
trajectory passes through a neighborhood of the index k saddle. Whether or not a given 
coordinate qi changes sign depends on the value of the integral /». If qi changes sign in the 
vicinity of the saddle then Jj > 0. If Jj = then q^ is zero or evolves to zero either in 
forward or negative time (depending on the initial condition). We consider the boundary of 
the region of crossing trajectories, characterised by Jj = in more detail below. 

More precisely, crossing and other trajectories are characterized as follows: 

(i) If Jj > for all i = 1, . . . , k, then the trajectory is a crossing trajectory. 

(ii) If Jj = and q^ = for all i = 1, k, then the trajectory is on the NHIM and is not 
a crossing trajectory. 

(iii) If Jj = for any i = 1, . . . ,k and qi ^ for any point on the trajectory then qi 
for all points on the trajectory (for finite time); the coordinate qi therefore does not 
change sign, and the trajectory is not a crossing trajectory. 

The signs of the integrals provide a rather coarse descriptor for crossing trajectories. In 
ref. [32] a symbolic description of saddle crossing trajectories for index k = 2 saddles was 
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introduced based on the sign change of the g^. This symbohc description distinguishes all 
quahtatively distinct classes of crossing trajectories. Here we recall the discussion of the 
symbolic description of crossing trajectories given in ref. [321 

Consider the index k = 2 case. The symbolic description of the behavior of a trajectory 
as it passes through a neighborhood of the index 2 saddle with respect to the coordinates 
(jk, k = 1,2, is expressed by the following four symbols, (/i/2;^i^2), where ii = ±, i2 = ±, 
/i = ±, /2 = ±. Here i^, k = 1,2, refers to the "initial" sign of as it enters the 
neighborhood of the index 2 saddle and fk, k = 1,2, refer to the "final" sign of q^, as it 

leaves the neighborhood of the index 2 saddle. For example, trajectories of type ( ; H — ) 

pass over the barrier from qi > to gi < 0, but remain on the side of the barrier with q2 < 0. 

Based on the number of distinct sequences of length four of + and — , there are 2^ = 16 
qualitatively distinct classes of trajectory. However, there are only four types of trajectories 
for which there is a change of sign of both coordinates qi and q2 as they pass through 

a neighborhood of the index 2 saddle, and these have symbolic descriptions (++; ), 

( — h; H — ), (H — ; — h) and ( ; ++). We will see that this symbohc description of crossing 

trajectories can be directly related to the geometry of the dividing surface. 

As noted previouslj^, codimension one surfaces separate the different types of trajectory. 
For index 2 saddles these are the codimension one invariant manifolds, given by qi = pi, 
^2 = P2 (i-e., /i = 0, J2 = 0). For example, the codimension one surface qi = pi forms the 
boundary between trajectories of type (++; H — ) and ( — h; H — ), and so on. 

An obvious generalization of this symbolic construction can be carried out for the case 
of index k saddles, k > 2, where there are 2"^^ possible classes of trajectory. 
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IV. THE DIVIDING SURFACE ASSOCIATED WITH INDEX k SADDLES 



An essential component of the analysis of reaction dynamics from the phase space per- 
spective is the construction of a dividing surface (DS) through which all reactant trajectories 
must pass, and which (locally) has the essential no-recrossing propertj^. In this section we 
discuss construction of a codimension one (within the energy surface) DS for an index k sad- 



dle using the equations of motion (2.8) associated with the normal form Hamiltonian of Sec. 
[TTJ This construction again demonstrates the utility of an analysis in terms of the 'physical' 
saddle coordinates gj. In addition to constructing the DS, using the normal form we are able 
to (locally) classify reactive trajectories, thereby obtaining a rigorous phase space character- 
ization of the subset of barrier crossing trajectories associated with concerted isomerization 
dynamics. 



A. Definition of dividing surface in the general case 

We define a measure of the distance from the origin (i.e., the saddle) in the configuration 
space of the saddle degrees-of-freedom as: 

k 

T^^\Y.^I (4-1) 

i=l 

The dividing surface we define contains the set of phase space points corresponding to the 
minimum distance from the origin attained by crossing trajectories. At any turning point 
in the variable T) we have 

k k 

'i^ = ^M^ = Y^ ^^liPi = 0- (4.2) 

1=1 i=l 

This equation defines a critical point for D, but it is in fact a minimum since 

k k 

i? = ^ A, (^,pi + qS^) = 5^ A- [ti + Pi) > 0, (4.3) 

i=l 1=1 

and where we have used J^Aj = for dynamics under the normal form Hamiltonian. Note 
that V = precisely when qi = pi = 0, i = 1, . . . , k, i.e., it is zero on the NHIM. 

The dividing surface at constant E is defined by the intersection of the following two 
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2n — 1 dimensional surfaces: 

5*1(^1,^1, . . . ,qk,Pk,qk+i,Pk+i, ■ ■ ■,qn,Pn) = H{Ii, ...,/„)- E = 0, (4.4a) 

k 

S2iqi,Pi, ■ ■ ■ ,qk,Pk,qk+i,Pk+i, ■ ■ ■,qn,Pn) = ^KqiPi = 0. (4.4b) 



i=l 



Points on the DS lie on crossing trajectories and must therefore satisfy the additional con- 
ditions: 

li>0, i = l,...,k. (4.5) 



The surface defined by (4.4) without the constraint in (4.5) includes trajectories which 
are not crossing trajectories and we will refer to it as the extended dividing surface. The 
boundary between the dividing surface and the rest of the extended dividing surface consists 
of phase points satisying the conditions Ji = 0, and/or I2 = and is codimension two, the 



intersection of (4.4) with the codimension one invariant manifolds Ji = 0, I2 = 0. 



It should be clear that the surface defined by (4.4) is codimension one, and is restricted 
to the energy surface H{Ii, In) = E hy construction. However, we need to prove that it 
has the properties of a DS, that is, 

(i) The vector field is transverse to the DS. 

The DS is codimension two in the full phase space (and codimension one restricted to 
the energy surface) and therefore has two vectors transverse to it at each point. 

The Hamiltonian vector field vh associated with the normal form Hamiltonian H is: 
vh = (AiPi,Aigi, . . . AfcPfe, Afcgfc,rifc+iPfc+i, -^k+m+i, ■ ■ ■ ,^nPn, -^nq-n) ■ (4.6) 

We will show that the Hamiltonian vector field Vh is transverse to the DS on the 
energy surface 5*1 = 0. 

The rate of change of 5*1 along vh is necessarily zero: 

dSi{vH) = {H,H} = 0, (4.7) 
where {■, ■} denotes the usual Poisson bracket. The rate of change of 5*2 along Vh is 

dS2{vH) = {S2,H} = j^S2 (4.8a) 



i=l 
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from (4.3). Clearly, this expression is greater than or equal to zero. It is zero precisely 



when Qj = pj = 0, j = 1, . . . , k, i.e. it is zero on the NHIM. 
(ii) All crossing trajectories pass through the dividing surface. 



It is evident that all trajectories of (2.8) that enter and leave a neighborhood of the 
origin necessarily achieve such a minimum distance from the origin. In particular, 
this is true for all trajectories satisfying Jj > 0, i = 1, . . . ,k. Therefore all crossing 
trajectories pass through the dividing surface. 

B. Index 1 Saddles 

It is instructive to see how our construction of a dividing surface for index k saddles 
reduces to the familiar case of index 1 saddleJ^Hl. in this case conditions (4.4) become 



Hih,...,Q = E (4.9a) 

qiPi = 0, (4.9b) 

where points on the dividing surface are also subject to the additional constraint: 

h = l {Pl - ql) > 0, (4.10) 

Now, in order to have qipi = and Ii = ^ {pi — qf) > we must have 

qi = 0, (4.11) 

which is the DS for index 1 saddles expressed in terms of the normal form coordinates q. 



C. Explicit parametrization of the DS for quadratic Hamiltonians 

For the case of a quadratic Hamiltonian it is possible to give an explicit parametrization 
of the dividing surface. For simplicity we consider the case of an index 2 saddle for a 
system with 2 DoF. For the two degree-of-freedom index 2 saddle the dividing surface is 2 
dimensional in the 3 dimensional energy surface, and visualization of the dividing surface 
is therefore possible. Such visualization provides insight into the crossing dynamics for an 
index 2 saddle, and we therefore discuss this topic in some detail. Moreover, we will see 
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that the dividing surface in the quadratic Hamiltonian case plays an important role in our 
algorithm for sampling the dividing surface associated with the full normal form for general 
(non-quadratic) Hamiltonians. 

The discussion given below is in the spirit of ref. [79l where the geometry associated with 
the quadratic Hamiltonian was used to to examine several different approaches to visualizing 
the phase space structures that govern reaction dynamics for index 1 saddles. A quadratic 
approximation to the Hamiltonian is expected to accurately capture the relevant geometrical 
features for energies 'close' to the energy of the saddle. 

We examine an n DoF quadratic Hamiltonian with an index 2 saddle as motivation for 
the parameterisation we will develop but the value of the index k in our formulas is kept 
general for use later. In the quadratic case, Aj = Aj = constant, Vli = uji = constant, and 



the Hamiltonian (2.2) reduces to: 



H{h,h, Is, ...,/„) = ^ XJi + ^J-^i (4-12a) 

i=l j=k+l 

= E iP^ - ^1) + E + (4-12b) 

i=l j=k+l 

which defines the energy surface: 

H{h,h,h,...,Q=E, (4.13) 

where we will consider the case E > 0. We define the energy Eg associated with the saddle 
DoF to be 

k n ^ 

Es = Y. iP" - ^1) = ^ - E fe' + P") ■ (4-14) 

i=l j=k+l 

(Aj and Ui constants) and rewrite this equation as follows: 

k k 



Es + J2l>^^^" = T.l^'^P" (4.15a) 

1=1 i=l 

k n -. 



2 ^2 

i=l j=k+l 



Two ellipsoids are defined in the phase space (g^ follows. For some 
parameter R> the condition 

k 

R = Y.\^^^'^ (4.16) 



i=l 
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defines a (fc — l)-dimensional ellipsoidal surface in the q configuration space and 



Es + R 



k 



^"2 



(4.17) 



defines a (/c — l)-dimensional ellipsoidal surface in the p momentum space. 

For the quadratic Hamiltonian, the condition for a trajectory to achieve a minimum 



distance from the origin in the q configuration space is (cf. (4.2)) 



^ = J] g^gi = J] \qiPi = (4.18) 

i=l i=l 

The case of an index 2 saddle with no center degrees-of-freedom corresponds to a two 
degree-of-freedom quadratic Hamiltonian system having an equilibrium point at the origin, 
where the matrix associated with the Hamiltonian vector field has two pairs of purely real 
eigenvalues, ±Ai, ±A2, with Ai,A2 > 0. The two dimensional dividing surface in the three 
dimensional energy surface given hj H = Es = E is then given parametrically by: 

sin 6* cos 6*' 



(91,92) 



'2R 



ip,,P2) = ±V2iE + R)(^ 



A2, 
cos 9 



sin 6 



(4.19a) 
(4.19b) 



Ai V A2 

where < i? < 00, < 6 < 2tt. It is straightforward to check that points on surface 



(4.19) satisfy conditions (4.16), (4.17), and (4.18), and that the Hamiltonian vector field is 



transverse to this surface. However, as we have not yet imposed condition (4.5), the extended 



dividing surface (4.19) contains phase points in addition to those on crossing trajectories. 



Projections of the full surface (4.19) into various 3 dimensional subspaces of phase space 



are shown in Figure [2| without the constraints on R and 6 implied by the conditions (4.5). 
Consequently, not all phase points on this surface lie on crossing trajectories, but the true 
DS is embedded in it. The parameter values chosen are Ai = 1, A2 = a/3, E = 1.0, and we 



show both signs of the square root in (4.19), with < -R < 1. 



Figures [2^ and [2]b show pi and p2, respectively, as functions of qi and ^2- Note that 
Pi, P2 are apparently continuous as pi, p2 pass through zero. This is however deceptive. 



If we examine equation (4.19) we see that for any configuration (gi,g2), the two parts of 
the surface given by the positive and negative square roots in the second equation must be 



distinct, since from equation (4.17) even when R is zero the values of pi, p2 cannot both be 



zero. Hence, although in Figure [2^ the values of pi are continuous on the axis ^2 = 0, pi = 
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the values of p2 are not continuous as can be seen from Figure |2)d, so that the two surface 
components are disjoint. 



Also note that, following the discussion in Section III, there are codimension one surfaces 
separating the different types of trajectory. These are the four codimension one invariant 
manifolds defined by particular values of the integrals and given by qi = pi, qi = —pi, 
q2 = p2, q2 = —p2 (i-e., h = 0, I2 = 0). For example, the codimension one surface qi = pi 
forms the boundary between trajectories of type (++; H — ) and ( — h; H — ), and so on. 

Figure |3] shows only the parts of the dividing surface satisfying the constraints given by 



equation (4.5). In other words, we only plot points on (4.19) for which Ii > 0, I2 > 0. 

In order to obtain explicit constraints we substitute the expressions for the phase points 
in terms of R and 6 into the integrals to obtain 



Ai Ai 



{E + 2R) _ 2 



and 



Ai 



cos' 9 



2R 

a7 



(4.20a) 
(4.20b) 

(4.20c) 



^sin2^^+^(sin2^^-cos2^^) 
A2 A2 



{E + 2R) 2 



A, 



sin'^ 



A^' 



(4.21a) 
(4.21b) 

(4.21c) 



Imposing the conditions Ji > and /2 > we then obtain the constraint 



R 



E + 2R 



< sin^ 9 < 



E + R 
E + 2R 



(4.22) 



Any point {R, 9) on the dividing surface satisfying condition (4.22) lies on a crossing trajec- 
tory. 

Figures [3^-[3jl show four almost disjoint, apparently continuous components of the surface 
which meet only at points pi = 0, p2 = where i? = (gi = g2 = 0) and pi, p2 satisfy 



equation (4.17). 



It is apparent from equation (4.19) that, for points on the DS with < < | (quadrant 
^1 > 0, g2 > 0), for the positive square root we have pi > 0, p2 < and CC trajectories 
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of the type (H — ; — h), while for the negative square root we have pi < 0, p2 > and CC 
trajectories of the type ( — h; H — ). 

Similarly for tt < < =y with points in the quadrant qi < 0, q2 < 0, for the positive square 
root we have pi < 0, p2 > and trajectories of type ( — h; H — ), while for the negative square 
root we have pi > 0, p2 < and points corresponding to CC trajectories of type (H — ; — h). 
There are therefore four disjoint pieces of the dividing surface, two of type (H — ; — h) and 
two of type ( — h; H — )• 

For the angle ranges | < ^ < vr and ^ < 6 < 27i there are similarly four disjoint pieces 
of the dividing surface: two of the type (++; ) and two of type ( ; ++). 

Figures and |4|d show pi{qi,q2) and ^2(^1, ^2), respectively, for the portion of the DS 
corresponding to CC trajectories of type (++; ). 



D. Sampling the dividing surface 

Sampling the DS for an index 2 saddle in the NF coordinate set and then using the sample 
points as initial conditions enables us to obtain directly concerted crossing trajectories of 
particular dynamical (symbolic) type. (See Sec. fV}) 

Although we have obtained a parametrization of the DS for the index 2 saddle for the 
quadratic Hamiltonian case, we would like to be able to sample the DS for general Hamil- 
tonians. Sampling the DS in the quadratic case is in principle straightforward: we need to 



sample points parametrized by {R, 6) in eq. (4.19) and use the NF coordinate transformation 



to convert these points to the original set of physical coordinates. The resulting phase points 
can be used as initial conditions for trajectory integration. 

For the case of a general Hamiltonian, there are however several complications: 



i) The frequencies Ai, A2 in (4.19) are no longer constant, but are in general functions 



of the phase space coordinates. 



ii) The energy E in (4.19) is in general a nonlinear function of the action variables. 



(iii) The truncated NF coordinate transformations and the associated Hamiltonian are 
approximate, and become less accurate the farther we move from the saddle point. 

(iv) In addition to sampling the saddle coordinates, for n > 3 DoF it is necessary to 
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sample the center planes, essentially partitioning the total energy between the center 
and saddle degrees of freedom but subject to the constraint H = E. 

(v) To use the sampled points as initial conditions for trajectories in microcanonical (con- 
stant E) simulations we require them to lie on the energy surface of the full Hamilto- 
nian. 

(vi) Ultimately, one would like to sample the DS uniformly with respect to some prescribed 
probability density. 

For the full nonlinear normal form one of the difficulties is that we need to sample the 
actions li, i = 1, . . . , k, and Ij, j = k + 1, . . . , n, subject to the nonlinear energy constraint 



(4.13) 



First, consider the problem of sampling phase points in the center planes. The phase 
space of the center DoF can be sampled uniformly, either in a rectangular grid in physical 
coordinates {qj,pj) or in action-angle variables but we need to restrict the range of 



these variables so that the constraints implied by the condition (4.13) on the total energy 
are obeyed. 

There are two possibilities. 

First, with zero energy in the saddle planes (Jj = 0, i = 1,2), we could calculate the 



allowed range for each phase space coordinate pair {qj,Pj) or each action Ij from (4.13). 
These values are however non-linearly interdependent and the calculation generally involves 
an iterative numerical procedure such as Newton's method. 

Alternatively, we can sample the center mode phase space over a sufficiently large but 
fixed region and then simply accept points for which H{Ii, . . . ,In) < E with Jj = for 
i = 1, . . . ,k. 

For each set of center DoF variables so obtained we need to sample points in the saddle 
planes such that H{Ii,. . . ,In) = E. For k = 1 (index 1) there are no free saddle plane 
parameters. For index k = 2 we sample values for parameters 6 and R. It is then necessary 



to solve (4.13) and (4.18) numerically for the values of {qi,Pi) 



For the general (non-quadratic) index 2 saddle we see that a point with given {R, 9) in 



(4.19) lies (by construction) on the surface 5*1, but does not in general lie on 5*2 (although 



it might approximately do so). Such a phase point satisfies (4.2), but (4.19) now does not 
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include the nonlinear (higher order) terms in actions /i . . . , /„ in an expansion of H . We 
define the energy correction via the relation 

k n 

Hih, h, /3, ...,/„) = 5^ AJ, + + (4.23a) 



i=l 
k 



j=k+l 



(4.23b) 



i=l 



j=k+l 



where AE includes the higher order terms from (4.19) but will also include NF errors which 
are in general functions of qi,pi. 



We now postulate a generalised parameterisation of (4.19) appropriate for an index 2 
saddle with n — 2 bath modes: 

sin 5 cos 6' 



(91,^2) 

(Pl,P2) 



/AT 



±V2(^. + i?) 



A2 

cos 6 — sin 6 



(4.24a) 
(4.24b) 



where R > 0, < 6 < 2tt. Note that the Aj are in general no longer constants so that 



the relations (4.24) are implicit equations for q = (gi, . . . , g^) and p = {pi, . . . ,Pk)- It is a 



simple matter to check that phase points parametrized by (4.24) satisfy generalised forms 
of (4.16), (4.17), and (4.18) (with the A replaced by A) if we note that ^^=1 Aj Jj formally 



defines the saddle energy Eg, with 

k 



i=l 



E ^^^^ 

j=k+i 



AE. 



(4.25) 



Such phase points satisfy both equations (4.2) and (4.4) and so lie on Si and S'2. 



For given parameters R,d the phase space coordinates {qi,Pi) must be calculated in 
combination with the Ai,Es, which are non-linear functions of the {qi,Pi)- While this is in 
general a difficult task, for the calculations reported here, the required values can be found 
by an iterative procedure based on the quadratic approximation to the Hamiltonian. 



Extension of eq. (4.19) to the index k saddle case with n — k bath modes is possible but 
cumbersome. Although we do not go into details here, a similar technique can be applied. 

Once we have a parametrization of the DS, we can use it to choose points on the surface 
and integrate associated trajectories forward and backward in time. Quantitative calculation 
of, for example, fiuxes requires sampling the DS according to a prescribed density or properly 
weighting the samples. 



21 



V. INDEX-2 SADDLES: MODEL POTENTIALS 



In this section we consider isomerization dynamics in a 2 DoF model potential exhibiting 
an index 2 saddle. Using the definition of the dividing surface in the vicinity of the index 
2 saddle discussed above, we are able to sample phase points on crossing trajectories via 
the normal form and the associated symbolic code. It is found that the trajectories defined 
in this way have the dynamical attributes one would intuitively associate with trajectories 
following a concerted isomerization mechanism; in this section we therefore refer to such 
trajectories as concerted crossing (CC) trajectories. 



A. 2 DoF 4 well model potential 

We consider a non-separable 2 DoF 4-well potential of the form 

v{Qi, Q2) = -aQl + Qt-Ql + Qt + PQlQl (5.1) 

In our numerical computations we use parameter values a = 2, /3 = 0.4. Contours of the 
potential for these parameters are shown in Fig. |5| Note that in this figure the horizontal 
axis is Q2 and the vertical axis is Qi. 
The associated Hamiltonian is: 

H=^ + ^ + viQ,,Q2), (5.2) 

with equations of motion: 

Qi = gp^= (5-3a) 

Q^=gp^= ^2, (5.3b) 

P, = - = 2aQ, - AQl - 2/3Q^Ql, (5.3c) 
oQi 

f) M 

^2 = = 2O2 - ^Ql - 2/3QIQ2. (5.3d) 

dQ2 

For the range of parameter values of interest, the potential v has the following set of 
critical points: 

(i) Minima (4) 



W-^^)=I±W^,±\/^K v = - (5.4) 
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(ii) Index-1 saddle (2) 



(iii) Index-1 saddle (2) 



(iv) Index-2 saddle (1) 



(Qi,Q2)= (±^,0J, v = -- (5.5) 



(gi,Q2)= (0,±--^^, v = -\ (5.6) 



(Qi,Q2) = (0,0), t; = 0. (5.7) 



The basic problem of interest associated with a potential of the form (5.1) concerns 
the nature of the isomerization dynamics (i.e., well-to- well transformations). In particular, 
we wish to distinguish in a dynamically rigorous and useful way between 'concerted' and 
'sequential' isomerization processes. 

Taking a > (5/2 for definiteness, we have the following rough classification of dynamics 
as a function of energy: 

• Confined regime: 

Trajectories are trapped in the vicinity of one of 4 possible minima, and no isomeriza- 
tion is possible. 

• Restricted isomerization: 

tt^ 1 , , 

-T<^^4 ^'-'^ 

Trajectories can pass between pairs of wells connected by the lowest energy index 1 
saddles. Passage between wells connected by higher energy saddles is not possible. 

In this case the standard phase space picture can be applied to analyze isomerizations, 
with reactant regions identified in the usual way. There are 2 symmetry related NHIMs, 
and for each isomerization reaction either a 2-state (RRKM) or 3-state (Gray-Ric d^° * ^^ *) 
model can be used to describe the isomerization kinetics for each pair of wells. 

• Unrestricted isomerization: 

-^<E <0. (5.10) 

Trajectories have sufficient energy to pass from any well to any well, but do not have 
enough energy to reach the index-2 saddle. 
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In this regime, the only possible way for the system to pass from the lower left well to 
the upper right well, say, is via sequential isomerization routes that proceed through 
either of the intervening wells (top left or lower right). Concerted passage via hilltop 
crossing is not energetically feasible. 

In this case there are 4 relevant NHIMs. A standard RRKM model could presumably 
be be used to analyze the kinetics (flux over saddles), or a 5-state generalized Gray- 
Rice model can be appliecP^'^. In the latter case, there are 4 reactant regions (wells) 
with boundaries consising of broken separatrices, and a 5th region lying outside the 
reactant region. We do not pursue such an analysis here. 

• "Roaming/concerted crossing" regime: E >0. 

Trajectories can wander freely over a single connected region of configuration space 
that encompasses all four wells. In this regime, a direct concerted isomerization route 
exists that connects, for example, the lower left and upper right wells. In principle 
this route coexists with the sequential routes discussed above; the dynamical problem 
addressed here then concerns the possibility of making a rigorous distinction between 
trajectories associated with the concerted and sequential isomerization pathways. This 
classification will be made in phase space using the normal form computed in the 
vicinity of the index two saddle. 



B. Sampling trajectories on the DS 



We wish to sample trajectories having prescribed dynamical character (CC and non-CC) 



for the 2 DoF system with Hamiltonian (5.2) using the DS specified by (4.4) with crossing 



trajectories subject to the constraint (4.5). The CC trajectories constitute a dynamically 



well-defined subset of trajectories that are associated with concerted well-to- well transitions. 
Isomerizing trajectories that enter the vicinity of an index k saddle and which are not crossing 
(non-CC) trajectories are potentially associated with sequential well-to-well transformations. 



The DS is sampled as outlined in Section |IVD[ A regular grid in cartesian coordinates 
(^ = Rcos[9],ri = _Rsin[^^]) provides an associated set of {R,9) values; the quadratic fre- 



quencies Ai, A2 are substituted for the exact, nonlinear frequencies Ai, A2 in (4.19) as a first 



approximation to give a sample point in the space of NF coordinates. New values for Ai, A2 
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are then computed from (2.3), and simple iteration gives a point satisfying (4.19). 

In general, the phase point thus obtained is not on the energy surface, so we must solve 
numerically for the value of Eg for which the condition on the total energy is satisfied. This is 
done by appropriate scaling of the momenta p. Each calculation of the energy error involves 
a separate iteration to obtain the nonlinear frequencies A^; for the calculations reported here 
a single iteration proves sufficient to achieve convergence. 

A further difficulty arises as a consequence of the inherent inaccuracy of the truncated 



NF transformations. If we solve for phase points with fixed NF energy as defined by (4.13) 
and then use the NF coordinate transformation to map the point into the original phys- 
ical coordinates, the resulting point obtained no longer lies exactly on the energy surface 
defined by the original physical Hamiltonian; the associated error increases with distance 
from the saddle. We must therefore determine the value of Es yielding a NF point which, 
after transformation, lies on the energy surface defined by the original Hamiltonian (ex- 
pressed in the original physical coordinates). This procedure yields points which we use as 
initial conditions for integrating trajectories of the full Hamiltonian in the original physical 
coordinates. 

In Figures |6| [7| [s] we show trajectories in {Qi,Q2) space obtained by sampling points on 
the DS as described above. The spacing of the sampling grid in (^,77) is set at 0.01 with 
R < 0.1. Our sampling procedure leads to a sparsity of points on the DS near the origin in 
configuration space, {Qi, Q2) = (0, 0). Note also that 6 is undefined at i? = 0. 

The total energies for the 2 DoF system are set aX E = 0.01, 0.1, 0.5, respectively. As we 



do not impose condition (4.5), we sample the extended dividing surface and so obtain points 
in addition to those on concerted crossing trajectories. 

At each sample point a trajectory is integrated forwards and backwards in time, classified 
by symbolic code (as above) and plotted. Each of the Figures [6| [T] and |8] has 4 panels, 
where each panel shows trajectories associated with points on the DS having the following 

symbolic classifications: concerted crossing (classes (++; ) and (H — ; — +)), and non-CC 

trajectories (classes (H — ; ) and (++; — h)). 

The results shown in Figs [6]-[8] are noteworthy in several respects. First, it is clear that, 
for all three values of the energy considered, the symbolic classification of trajectories as 
obtained using the NF in the vicinity of the index-2 saddle corresponds precisely to the 
observed dynamical behavior of the numerically integrated trajectories. That is, the NF en- 
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ables us to sample a subset of trajectories of a given dynamical type, e.g., concerted crossing. 
Second, we note that, although the ensemble of CC trajectories for a given energy appears 
to consist of two disjoint pencils or bundles, the set of initial consitions is in fact connected 
(there is a CC trajectory passing through the origin). The form of the configuration space 
projections of the various trajectory classes is by no means inuitively obvious: at all three 
energies studied, CC trajectories tend to be concentrated away from the hilltop itself. The 
non-CC trajectories passing through the vicinity of the index 2 saddle tend to 'bounce' 
off the saddles as they pass between the wells, while CC trajectories appear to 'graze' the 
hilltop. 

As mentioned previously, the boundary between CC and non-CC trajectories on the 
extended dividing surface is composed of those points on the DS satisfying Ji = and/or 
I2 = 0. Projected into configuration space, the boundary consists of a number of lines 
emanating from the origin. In Fig. [9] we show the projection into configuration space of a 
segment of a boundary between CC and non-CC trajectories defined by the condition I2 = 
with < R < 0.1 and E = 0.01. Each phase point on the boundary is propagated forward 
and backward in time. It can be seen that trajectories associated with boundary points 
start in the lower left hand well, pass through the vicinity of the saddle and end up (at 
short times) tending towards the vertical axis, i.e., neither left (non-CC) nor right (CC). 
Again, the phase points on the boundary are obtained using the NF, yet show precisely the 
expected dynamical behavior when propagated numerically. 

Lastly, it is natural to consider the fraction of CC trajectories in a given ensemble at 
a prescribed energy. To provide an unambigous determination of this quantity, we must 
carefully define the relevant ensemble, such as we have done above in terms of the sampling 
procedure on the DS, and also specify the relevant weighting factor (measure) for trajectories. 
As the coordinates {R, 9) used to parametrize the DS are not canonical coordinates, a 
weighting factor enters that is a non-trivial function of {R,6). Rather than calculate the 
fraction of CC trajectories by sampling the DS using this weight function, in the next 
subsection we consider a different sampling procedure where trajectories are initiated on the 
plane Qi = with Pi > 0. In this case the associated density at constant energy is just the 
natural measure (area) in the {Q2, P2) plane. 
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C. Trajectory studies of isomerization dynamics 



We now study the isomerization dynamics in the concerted crossing regime using an 
approach complementary to that of the previous subsection. These computations provide 
additional insight into the way in which the presence of the index-2 saddle affects trajectories 
in the neighborhood of the saddle, and further confirm the accuracy of the NF in the vicinity 
of the index 2 saddle. 

Initial conditions are chosen as follows. In the configuration space {Qi,Q2)j we define a 
regular grid of points along the line (Qi = in an interval symmetric about Q2 = 0. At each 
point {Qi,Q2) along this line, a number of momentum pairs (Pi,P2) are chosen such that 
the phase space points {Qi,Q2, Pi, A) he on the fixed energy surface, H = E. We consider 
the same 3 energies as in the previous subsection, H{Qi, Q2, -Pi, A) = 0.01, 0.1 and 0.5, and 
we set dQi/dt > 0, so that trajectories start out moving from the lower to the upper half 
of the potential (the potential is symmetric by construction). 

A symbolic code can be assigned to each initial condition in either of two ways. Tra- 
jectories can be integrated forward and backward in time until the first well is reached. In 
practice this means that we integrate trajectories until a turning point is reached in either of 
the coordinates (5i,2- In forward time the trajectory can enter one of 2 possible wells, and in 
backward time the trajectory can enter 2 possible wells, so that there are 4 qualitatively dif- 
ferent type of trajectories. Using our symbolic classificatioEp', the trajectory (/i, /2; ^1, "^2), 
where = ±, = ±, = 1,2, denotes the trajectory that passes from well (ii,«2) in the 
immediate past to the well (/i,/2) in the immediate future. Each initial condition on the 
horizontal line can therefore be labelled according to its symbolic description (/i, /2; ^1, ^2), 
i.e., the first well visited in the past and future, as determined by the exact trajectory dy- 
namics. We can also use the NF to classify phase points in the Qi = plane, and this 
symbolic classification can be compared with the trajectory results. 

Initial conditions at constant energy with Qi = 0, Pi > are uniquely specified by the 
values of the phase space variables ((52,-^2)- We therefore assign the symbolic trajectory 
code (4 possibilities) as a function of coordinates ((52,-^2)- Initial conditions are sampled 
uniformly in the ((52,-^2) plane. 



Our results are presented in Fig. 10 in which we plot the fraction F of crossing trajectories 
as a function of coordinate Q2- Results are shown for three energies E = 0.01, 0.1 and 
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0.5 across the interval —0.5 < Q2 < 0.5. This range of the coordinate Q2 corresponds 
approximately to the constraint < < 0.1 employed when sampling the DS as described in 
the previous subsection. On the same graph we show the corresponding trajectory fractions 
determined by taking each initial condition, converting it to NF coordinates, and using the 



NF to predict the symbolic trajectory code. In Figure [10] the sample spacing for trajectory 
initialisation was /S.Q2 = 0.01, while for the faster normal form sampling we used a 5 times 
finer spacing, AQ2 = 0.002. Use of the NF results in a smoother curve. If the same 
sampling spacing is used as for the trajectory calculations, the NF results are essentially 
identical to those obtained by trajectory integration. Discrepencies between the trajectory 
and NF results are most pronounced near the boundaries of the interval, and this is not 
unexpected as the accuracy of the NF presumably deteriorates as one moves away from the 
saddle. What is perhaps striking about our results is the size of the region of the Q2 axis 
over which the NF predictions do accurately match the trajectory results. 



D. Saddle crossing in the presence of bath modes 

Finally, we briefly explore the nature of the saddle crossing dynamics in the presence of 
additional degrees of freedom (bath modes). We consider a 4 DoF model in which 2 bath 
modes are added to the 4-well 2 DoF system studied above, and bilinear coupling terms are 
introduced between saddle and bath modes. 



1. System-bath Hamiltonian 



We consider a 4 DoF system with an index 2 saddle point. The relevant system-bath 
Hamiltonian Hgh is obtained by adding two 'bath' modes, coordinates {x,y), to the 2 DoF 



system with Hamiltonian given by eq. (|5.2|): 



sb 



' ' V{Qi,Q2) + 



2 , / ciQi 

Px + ^xX 



+ 



2 , / C2(52 
Py + ^yV 



(5.11) 



where mode x is coupled to system coordinate Qi and mode y is coupled to system coordinate 
Q2. 
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Hamilton's equations are given by: 



dH 

Qi = ww = A, (5.12a) 
oPi 

Q^ = QP^= ^2, (5.12b) 
dH 

2;=^=Px, (5.12c) 
dH 

y = ^=Py, (5.12d) 



dp. 



y 



OH -3 ^nr^ Ci / CiQi 



P, = - = 2aQi - 4Qi - 2/3Q,Qi + ^ cu^x - ^ , (5.12e) 

oQi \ U)^ ' 



^2 = = 2Q2 - ^Ql - 2PQ\Q, + ^ Lyy - , (5.12f) 

0Q2 \ ' ^y J 

Px- = --^ = -uJi I — j , (5.12g) 

= = -^y ^yV • (5.12h) 



dy " ^y / 

For our numerical calculations we use bath parameters cUx = 1.0, cOy = \/2 and set Ci = C2. 

Sampling trajectories on the extended DS for 4 DoF 

We now consider the effect of additional degrees of freedom on the dynamics, in particular, 
our ability to sample concerted crossing trajectories on the (extended) DS for the 4 DoF 



system-bath Hamiltonian, eq. (5.11). 



The energy shell for the 4 DoF system-bath model is 7 dimensional, while the DS is 
6 dimensional. Fully sampling the DS according to the procedures described above is a 
numerically intensive task. For purposes of illustration, we demonstrate that computation 
of the NF for the 4 DoF model allows us to sample CC trajectories on the DS, and that 
bundles of such trajectories, when projected into the {Qi,Q2) subspace, appear simply as 
'fattened' versions of the 2DoF saddle trajectories described above with an appropriate value 
of the saddle energy. This inherent apparent simplicity of the system-bath dynamics is a 
consequence of the integrability of the NF in the vicinity of the saddle. 

In order to sample the DS, we first set the saddle energy Eg = (no excitation in center 
modes) and sample phase space coordinates in the saddle planes; for each saddle plane 
point we sample phase points in the center DoF setting x = ±0.5, y = ±0.5, = ±0.5, 
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and Py = ±0.5, and scale the center DoF coordinates and momenta to obtain a fixed value 
of the total energy i?sb = 0.5. 



Figure 11 shows {Qi,Q2) projections of trajectories obtained using the sampling proce- 
dure just described with saddle-bath coupling parameter ci = C2 = 0.1. 

Each initial condition in the saddle planes is therefore associated with a 'bundle' of 
trajectories for the 4 DoF full system. 



Figure [12] shows analogous results for larger system-bath coupling parameters Ci = C2 = 
0.5. For the larger coupling we obtain a thicker bundle associated with the fiducial initial 
condition in the saddle planes. 
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VI. SUMMARY AND CONCLUSIONS 



In this paper we have extended our earher analysis of the phase space structure in the 
vicinity of an equihbrium point associated with an index k saddled^. We have shown that 
Poincare-Birkhoff normal form theory provides a constructive procedure for obtaining an 
integrable approximation to the full Hamiltonian in the vicinity of the equilibrium, provided 
a generic non-resonance condition is satisfied independently for both the real eigenvalues 
and the complex eigenvalues of the matrix associated with the linearization of Hamilton's 
equations about the index k saddle. As a consequence there are k independent integrals 
associated with the saddle degrees-of-freedom. These integrals provide a precise tool for 
classifying trajectories that pass through a neighborhood of the saddle. In particular, they 
provide a symbolic classification of the trajectories into 2"^^ distinct types of trajectory that 
pass through a neighborhood of the index k saddle. 

The normal form also provides an algorithm for constructing a dividing surface, i.e., a 
co-dimension one surface (restricted to the energy surface) through which all trajectories 
in a neighborhood of the index k saddle must pass (with the exception of a set of zero 
measure). We provide a parametrization of this dividing surface which, when using the 
integrals associated with the saddle degrees-of-freedom, can be sampled in such a way that 
we can choose initial conditions corresponding to any particular type of trajectory described 
by the symbolic classification. 

We illustrated our analytical and computational techniques by analyzing a problem that 
brings to light a fundamental mechanistic role played by index two saddles in chemical 
dynamics. Namely, we consider isomerization on a potential energy surface with multi- 
ple symmetry equivalent minima. In the two degree-of-freedom example we computed the 
normal form and the dividing surface and showed that the different classes of reactive trajec- 
tories in the vicinity of the index two saddle (for three different energies) could be computed 
by our sampling routine. Our procedure enables a rigorous definition of concerted crossing 
trajectories to be given in terms of local phase space structure. We then considered a sim- 
plified system-bath model (one harmonic oscillator mode for each degree-of-freedom), and 
showed that our approach could be applied to this four degree-of-freedom system. 
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Figure captions 



FIG. 1: Schematic representation of sequential versus concerted isomerizing pathways in a 
model 4-well potential, (a) Sequential, (b) Concerted. 

FIG. 2: Projections of the index 2 saddle phase space extended dividing surface, parameter 
values Ai = 1, A2 = v^, E = 1.0, with < i? < 1. (a) Coordinates (gi, q2, pi)- (b) 
Coordinates (gi, q2, ^2)- (c) Coordinates {pi, p2, qi). (d) Coordinates {pi, p2, (12)- 
FIG. 3: Projections of the index 2 saddle phase space dividing surface with the constraints 



(4.5), parameter values Ai = 1, A2 = -\/3, E = 1.0. (a) Coordinates (gi, q2, Pi). (b) 
Coordinates {qi, q2, p2)- (c) Coordinates {pi, p2, qi). (d) Coordinates {pi, p2, (12) ■ 
FIG. 4: Projections of the portion of the index 2 saddle phase space dividing surface as- 
sociated with trajectories of symbolic type + H . Parameter values Ai = 1, A2 = v^, 

E = 1.0. (a) Coordinates {qi, q2, pi)- (b) Coordinates {qi, q2, p2)- 



FIG. 5: Contour plot of nonseparable 2 DoF potential v{Qi, Q2), eq. (5.1 ). Parameter values 
a = 2, 13 = 0.4. 

FIG. 6: Trajectories initiated on the DS and propagated forwards and backwards in time. 

Saddle energy is 0.01. (a) Concerted crossing trajectories (++; ). (b) Concerted crossing 

trajectories (H — ; — h). (c) Non-CC trajectories (H — ; ), I2 < 0. (d) Non-CC trajectories 

(++;-+), /2<0. 

FIG. 7: Trajectories initiated on the DS and propagated forwards and backwards in time. 

Saddle energy is 0.1. (a) Concerted crossing trajectories (++; ). (b) Concerted crossing 

trajectories (H — ; — h). (c) Non-CC trajectories (H — ; ), I2 < 0. (d) Non-CC trajectories 

(++;-+), J2<0. 

FIG. 8: Trajectories initiated on the DS and propagated forwards and backwards in time. 

Saddle energy is 0.5. (a) Concerted crossing trajectories (++; ). (b) Concerted crossing 

trajectories (H — ; — h). (c) Non-CC trajectories (H — ; ), I2 < 0. (d) Non-CC trajectories 

(++;-+), /2<0. 

FIG. 9: Segment of the boundary between CC and non-CC trajectories defined by the 
condition I2 = 0. Crosses indicate configuration space projections of boundary points with 
< -R < 0.1 and E = 0.01. Each phase point on the boundary is propagated forward and 
backward in time. 
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FIG. 10: Fraction F of crossing trajectories for initial conditions at constant energy in the 
Qi = plane with —0.5 < Q2 < +0.5. NF predictions (green hne) are show together with 
results obtained from integration of trajectories (red line), (a) E — 0.01, (b) E — 0.1, (c) 
E^O.5. 

FIG. 11: Trajectories initiated on the DS and propagated forwards and backwards in time. 
Total energy E — 0.1. Bath coupling parameter ci — C2 — 0.1. (a) Goncerted crossing tra- 
jectories (++; ). (b) Goncerted crossing trajectories (H — ; — h). (c) Non-GG trajectories 

(H — ; ), h < 0. (d) Non-CG trajectories (++; — h), h < 0. 

FIG. 12: Trajectories initiated on the DS and propagated forwards and backwards in time. 
Total energy E — 0.1. Bath coupling parameter ci — C2 — 0.5. (a) Goncerted crossing tra- 
jectories (++; ). (b) Goncerted crossing trajectories (H — ; — h). (c) Non-GG trajectories 

(H — ; ), /2 < 0. (d) Non-GG trajectories (++; — h), h < 0. 
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FIGURE 6 
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FIGURE 9 
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FIGURE 11 



48 




-1.5 -1 -0.5 0.5 1 1.5 -1.5 -1 -0.5 0.5 1 1.5 
(a) Q2 (b) Q2 




-1.5 -1 -0.5 0.5 1 1.5 -1.5 -1 -0.5 0.5 1 1.5 
(c) Q2 (d) O2 



FIGURE 12 



49 



